Riassunto Analisi Predittiva

Introduzione

Questo documento fornisce un’ampia panoramica sul corso di Analisi Predittiva di Ca’ Foscari (CT0429) dell’anno 2023/2024. Nel corso del documento, esploreremo vari argomenti riguardanti l’Analisi Predittiva, compresi Anova, Residui, Formule, Correlazione tra variabili, Predict, l’utilizzo di più modelli, e come selezionare le variabili all’interno di un modello.

L’obiettivo di questo documento è guidarti attraverso i concetti chiave e le pratiche nell’Analisi Predittiva, con un focus sull’utilizzo del linguaggio di programmazione R per applicare queste tecniche. Spero che questo documento ti aiuti a comprendere meglio questa materia di studio al fine di superare l’esame.

Indice

• Regressione Lineare
  • Analisi dei residui
  • Scomposizione Somma Di Quadrati
  • Predizione Ottimale
  • Covarianza e Correlazione Empirica
• Regressione Lineare Multipla
  • Approcio Matriciale
  • Interpretazione Geometrica
  • Distribuzione F & Anova Table
  • Modelli Nidificati
  • Selezione delle variabili
  • Predizioni Categoriche
  • Controllo del Modello
  • Transformazioni
  • Multicollinearità
  • Punti Influenti
• Modelli Lineari Generalizzati
  • Devianza
  • Residui

Regressione Lineare

La Regressione Lineare (LR) è una tecnica statistica che viene utilizzata per studiare la relazione tra due o più variabili quantitative. Il modello di regressione lineare assume che la variabile dipendente, che si desidera predire, sia una funzione lineare delle variabili indipendenti.

Il modello di regressione lineare fa le seguenti assunzioni:

• Linearità: La relazione tra la variabile dipendente e le variabili indipendenti è lineare.
• Normalità: I residui sono distribuiti normalmente.
• Indipendenza: I residui sono indipendenti tra loro.
• Eteroschedasticità: La varianza dei residui è costante.

L’equazione del modello di regressione lineare semplice è la seguente:

\[ y = a + bx \]

dove:

• y è la variabile dipendente
• a è l’intercetta
• b è il coefficiente angolare
• x è la variabile indipendente

L’intercetta rappresenta il valore medio di y quando x è uguale a 0. Il coefficiente angolare rappresenta la variazione di y per ogni unità di variazione di x.

Analisi dei residui

Serve a verificare se il modello soddisfa le assunzioni della regressione e per identificare eventuali pattern o problemi nei dati.

Esempio 1:

# Generiamo dati casuali con residui normalmente distribuiti
set.seed(123)
x <- 1:100
y <- 2 * x + rnorm(100)

# Adattiamo un modello di regressione
model <- lm(y ~ x)

# Effettuiamo l'analisi dei residui
Residuals <- residuals(model)

# Creiamo un grafico dei residui
plot(x, Residuals, main = "Distribuzione Normale dei Residui",xlab = "X", ylab = "Residui")
abline(h = 0, col = "red")

In questo caso, i residui seguono una distribuzione normale, il che è un risultato ideale per un modello di regressione lineare. Quando diciamo che i residui hanno una distribuzione normale, significa che i residui seguono una distribuzione a forma di campana, con una media di zero e una varianza costante. Questa è un’importante assunzione nei modelli di regressione lineare, in quanto indica che gli errori casuali nel modello sono distribuiti in modo simmetrico intorno a zero e non mostrano alcun tipo di tendenza sistemica. Se questa assunzione è soddisfatta, i test di significatività dei coefficienti del modello e le stime di intervallo di confidenza saranno affidabili.

Esempio 2:

# Generiamo dati casuali con residui che seguono una distribuzione a U
set.seed(456)
x <- 1:100
y <- 2 * x^2 + rnorm(100)
y[50:60] <- y[50:60] + 10  # Introduciamo un effetto a U nei dati

# Adattiamo un modello di regressione
model <- lm(y ~ x)

# Effettuiamo l'analisi dei residui
residuals <- residuals(model)

# Creiamo un grafico dei residui
plot(x, residuals, main = "Distribuzione a U dei Residui", xlab = "X", ylab = "Residui")
abline(h = 0, col = "red")

In questo caso, i residui mostrano un effetto a U, indicando una violazione dell’assunzione di omoschedasticità(i residui non mostrano un aumento o una diminuzione sistematica nella dispersione al variare dei valori delle variabili indipendenti).

Esempio 3:

# Generiamo dati casuali con outlier nei residui
set.seed(789)
x <- 1:100
y <- 2 * x + rnorm(100)
y[c(20, 85)] <- y[c(20, 85)] + 20  # Aggiungiamo outlier nei dati

# Adattiamo un modello di regressione
model <- lm(y ~ x)

# Effettuiamo l'analisi dei residui
residuals <- residuals(model)

# Creiamo un grafico dei residui
plot(x, residuals, main = "Presenza di Outlier nei Residui", ylab = "Residui", xlab = "X")
abline(h = 0, col = "red")

In questo caso, i residui mostrano la presenza di outlier evidenti, che possono influenzare in modo significativo la stima dei coefficienti del modello. Gli outlier possono comportare problemi nei modelli statistici, specialmente nei modelli di regressione, perché possono influenzare notevolmente i risultati. Ad esempio, possono influenzare la stima dei coefficienti del modello e rendere il modello meno affidabile.

[Torna all’ Indice]

Scomposizione Somma Di Quadrati

La “decomposition of sum of squares” è un concetto fondamentale nell’analisi della varianza (ANOVA) e nella regressione statistica. Questa tecnica aiuta a scomporre la varianza totale osservata in un insieme di dati in diverse componenti, consentendo di comprendere quanto della varianza può essere attribuito a vari fattori o errori residui. La formula chiave in questo contesto è:

\[ Varianza Totale = Varianza Spiegata + Varianza Residua \]

Dove:

• Varianza Totale è la varianza complessiva dei dati, cioè quanto i dati variano in generale.

• Varianza Spiegata rappresenta la varianza dovuta al modello o ai fattori esaminati (spiegati dalla variabile indipendente nel contesto della regressione).

• Varianza Residua è la varianza non spiegata dal modello o dai fattori ed è associata all’errore residuo, ovvero la differenza tra i valori osservati e quelli previsti dal modello.

Nel contesto della regressione, puoi rappresentare la decomposizione della somma dei quadrati come segue:

\[ SST = SSR + SSE \]

• SST (Sum of Squares Total) rappresenta la somma dei quadrati totale ed è la varianza dei dati osservati rispetto alla loro media.

• SSR (Sum of Squares Regression) rappresenta la varianza spiegata dal modello o dalla variabile indipendente.

• SSE (Sum of Squares Error) rappresenta la varianza residua, ossia la varianza non spiegata dal modello.

Per valutare l’efficienza del tuo modello di regressione, dovresti guardare la proporzione di questa varianza spiegata dal tuo modello (SSR). In generale, vuoi massimizzare la proporzione spiegata e minimizzare la proporzione non spiegata (SSE). Pertanto, punti a minimizzare SSE.

data <- data.frame(X = c(1, 2, 3, 4, 5), Y = c(3, 5, 6, 8, 10))
mean_Y <- mean(data$Y)
SST <- sum((data$Y - mean_Y)^2)

# Adatta il modello di regressione lineare
model <- lm(Y ~ X, data = data)

# Calcola la SSR
SSR <- sum((predict(model) - mean_Y)^2)

# Calcola la SSE
SSE <- sum(model$residuals^2)

R_squared <- SSR / SST

# Equivale a fare summary(model)
R_squared

[1] 0.989726

[Torna all’ Indice]

Predizione Ottimale

L’“Optimal Prediction” riguarda la determinazione di un modello predittivo che sia il migliore possibile in termini di accuratezza nel prevedere gli eventi futuri.

L’obiettivo principale è trovare il modello che massimizza la precisione delle previsioni, minimizzando l’errore di previsione. Ci sono vari metodi e tecniche per ottenere la predizione ottimale, a seconda del contesto e dei dati disponibili.

# Caricamento del dataset "cars"
data(cars)

# Visualizzazione delle prime righe del dataset
head(cars)

# Dividiamo il dataset in set di addestramento e set di test
set.seed(123)  # Impostiamo un seed per la riproducibilità
sample_indices <- sample(nrow(cars), nrow(cars) * 0.7)  # 70% dati di addestramento
train_data <- cars[sample_indices, ]
test_data <- cars[-sample_indices, ]

# Adattamento di un modello di regressione lineare
model <- lm(dist ~ speed, data = train_data)

# Predizioni
predictions <- predict(model, newdata = test_data)

# Valutazione delle prestazioni
summary(model)

Call:
lm(formula = dist ~ speed, data = train_data)

Residuals:
    Min      1Q  Median      3Q     Max 
-18.820  -8.798  -2.272   5.614  44.951 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) -22.0481     7.4169  -2.973  0.00548 ** 
speed         4.0457     0.4589   8.817 3.44e-10 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 13.08 on 33 degrees of freedom
Multiple R-squared:  0.702, Adjusted R-squared:  0.693 
F-statistic: 77.73 on 1 and 33 DF,  p-value: 3.435e-10

# Grafico dei risultati
library(ggplot2)
ggplot(data = test_data, aes(x = speed, y = dist)) +
  geom_point(color = "blue") +
  geom_smooth(method = "lm", se = FALSE, color = "red") +
  ggtitle("Predizione della Distanza di Arresto")

[Torna all’ Indice]

Covarianza e Correlazione Empirica

La covarianza e correlazione empirica sono misure statistiche utilizzate per quantificare la relazione tra due variabili in un insieme di dati osservati. Queste misure sono strettamente legate e sono spesso utilizzate per esaminare la relazione lineare tra due variabili.

Covarianza Empirica:

La covarianza empirica è una misura della tendenza di due variabili a variare insieme. Indica se le due variabili crescono o diminuiscono simultaneamente (covarianza positiva) o se una aumenta mentre l’altra diminuisce (covarianza negativa). La formula per calcolare la covarianza empirica tra due variabili X e Y in un set di dati è data da:

\[ Cov(X,Y) = \frac{1}{n-1} \sum_{i = 1}^n {(X_i - \overline{X})(Y_i - \overline{Y})} \]

Correlazione Empirica:

La correlazione empirica è una versione standardizzata della covarianza empirica e misura la forza e la direzione di una relazione lineare tra due variabili. La correlazione empirica è sempre compresa tra -1 e 1. La formula per calcolare la correlazione empirica tra due variabili X e Y è data da:

\[ Cor(X,Y) = \frac{Cov(X,Y)}{S_X \cdot S_Y} \]

Le misure di covarianza empirica e correlazione empirica sono utilizzate per esaminare la relazione tra variabili in un set di dati e sono particolarmente utili nell’analisi statistica e nell’apprendimento automatico per valutare le associazioni tra le variabili prima di costruire modelli predittivi. La correlazione empirica è più comunemente utilizzata perché fornisce una misura standardizzata della relazione tra variabili ed è meno influenzata dall’unità di misura.

Il coefficiente di Pearson è utile per vedere se due variabili hanno una correlazione lineare o meno. Questo perché non tutte le variabili correlate hanno una relazione lineare.

Esempio:

# Esempio dati casuali
set.seed(123)
x <- rnorm(100)  # Variabile x
y <- 2 * x + rnorm(100)  # Variabile y (correlata a x)

# Calcola la correlazione di Pearson
correlation <- cor(x, y)

# Stampa il valore di correlazione
cat("Correlazione di Pearson tra x e y:", correlation, "\n")

Correlazione di Pearson tra x e y: 0.8786993 

In questo esempio, stiamo generando dati casuali per le variabili x e y. La variabile y è costruita come una trasformazione lineare di x con un termine di errore aggiunto. Poi, utilizziamo la funzione cor() per calcolare la correlazione di Pearson tra x e y.

Un valore vicino a 1 indica una correlazione lineare positiva forte, un valore vicino a -1 indica una correlazione lineare negativa forte, mentre un valore vicino a 0 indica una scarsa correlazione lineare tra le due variabili.

[Torna all’ Indice]

Regressione Lineare Multipla

La “Multiple Linear Regression” (Regressione Lineare Multipla) è una tecnica di modellazione statistica utilizzata per analizzare la relazione tra una variabile dipendente (o target) e due o più variabili indipendenti (o predittive). Questa tecnica estende la semplice regressione lineare, che coinvolge solo una variabile indipendente, a un contesto in cui più variabili indipendenti sono coinvolte nel modello. La regressione lineare multipla è ampiamente utilizzata nell’analisi statistica e nell’apprendimento automatico per fare previsioni o comprendere le relazioni complesse tra variabili.

\[ Y = \beta_0 + \beta_1X_1 + ... + \beta_nX_n + \epsilon \] L’obiettivo principale è stimare i coefficienti β in modo che il modello si adatti meglio ai dati osservati. Questo viene fatto utilizzando metodi di stima, come il metodo dei minimi quadrati, che cerca di minimizzare la somma dei quadrati degli errori residui. Il modello di regressione viene valutato utilizzando metriche di valutazione delle prestazioni come l’errore quadratico medio (RMSE), il coefficiente di determinazione (R-squared) e altri. È importante eseguire test di significatività statistica per i coefficienti delle variabili indipendenti per determinare se esse contribuiscono significativamente al modello.

La regressione lineare multipla è basata su alcune assunzioni, tra cui l’indipendenza degli errori, l’omoschedasticità (varianza costante degli errori), la linearità della relazione e la normalità degli errori.

# Carica il dataset mtcars
data(mtcars)

# Visualizza le prime righe del dataset
head(mtcars)

# Adattamento del modello di regressione lineare multipla
model <- lm(mpg ~ wt + hp + qsec, data = mtcars)

# Visualizza un riepilogo del modello
summary(model)

Call:
lm(formula = mpg ~ wt + hp + qsec, data = mtcars)

Residuals:
    Min      1Q  Median      3Q     Max 
-3.8591 -1.6418 -0.4636  1.1940  5.6092 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) 27.61053    8.41993   3.279  0.00278 ** 
wt          -4.35880    0.75270  -5.791 3.22e-06 ***
hp          -0.01782    0.01498  -1.190  0.24418    
qsec         0.51083    0.43922   1.163  0.25463    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 2.578 on 28 degrees of freedom
Multiple R-squared:  0.8348,    Adjusted R-squared:  0.8171 
F-statistic: 47.15 on 3 and 28 DF,  p-value: 4.506e-11

# Effettua previsioni con il modello
predictions <- predict(model, newdata = mtcars)

# Valuta le prestazioni del modello
SSR <- sum((mtcars$mpg - predictions)^2)  # Somma dei quadrati residui
SST <- sum((mtcars$mpg - mean(mtcars$mpg))^2)  # Somma totale dei quadrati

R_squared <- 1 - (SSR / SST)

cat("R-squared (manual calculation):", R_squared, "\n")

R-squared (manual calculation): 0.8347678 

In sintesi, il modello di regressione lineare multipla suggerisce che il peso del veicolo influenza il consumo di carburante. La potenza del motore e il tempo di accelerazione non sono significative per la spiegazione di questo modello. Nonostante ciò il modello spiega l’83% dei della variazione nei consumi di carburante.

[Torna all’ Indice]

Approcio Matriciale

L’approccio matriciale alla regressione è una forma alternativa di rappresentazione e risoluzione dei modelli di regressione, inclusa la regressione lineare. Questo approccio utilizza notazioni matematiche e matrici per semplificare i calcoli e ottenere soluzioni più efficienti in problemi di regressione lineare.

Ci permette di passare da: \[ Y = \beta_0 + \beta_1X_1 + ... + \beta_nX_n + \epsilon \] a : \[ Y = \beta X + \epsilon \]

La soluzione matriciale per stimare i coefficienti β é: \[ \beta = ((X^T X)^{-1} X^T Y) \]

L’approccio matriciale semplifica la rappresentazione e la risoluzione dei modelli di regressione, specialmente quando si lavora con più variabili indipendenti. Inoltre, è utile per comprendere come eseguire calcoli di regressione in modo più efficiente utilizzando matrici e algebra lineare, specialmente in contesti di apprendimento automatico in cui le dimensioni dei dati possono essere elevate.

# Carica il dataset "swiss" (un dataset di dati demografici svizzeri)
data(swiss)

# Visualizza le prime righe del dataset
head(swiss)

# Crea la matrice delle variabili indipendenti
X <- as.matrix(swiss[, c("Examination", "Education")])

# Aggiungi una colonna di 1 per l'intercetta
X <- cbind(1, X)

# Crea il vettore delle variabili dipendenti
Y <- swiss$Fertility

# Calcola i coefficienti del modello utilizzando l'approccio matriciale
beta <- solve(t(X) %*% X) %*% t(X) %*% Y

# Visualizza i coefficienti del modello
print(beta)

                  [,1]
            85.2532753
Examination -0.5572183
Education   -0.5394570

# Effettua previsioni con il modello matriciale
predictions <- X %*% beta

# Valuta le prestazioni del modello
model <- lm(Y ~ Examination + Education, data = swiss)  # 0 indica di non calcolare l'intercetta
summary(model) 

Call:
lm(formula = Y ~ Examination + Education, data = swiss)

Residuals:
     Min       1Q   Median       3Q      Max 
-15.9935  -6.8894  -0.3621   7.1640  19.2634 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  85.2533     3.0855  27.630   <2e-16 ***
Examination  -0.5572     0.2319  -2.402   0.0206 *  
Education    -0.5395     0.1924  -2.803   0.0075 ** 
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 8.982 on 44 degrees of freedom
Multiple R-squared:  0.5055,    Adjusted R-squared:  0.483 
F-statistic: 22.49 on 2 and 44 DF,  p-value: 1.87e-07

• Esame (Examination): Un aumento nei punteggi di esame è associato a una diminuzione della fertilità nelle regioni svizzere. Questo suggerisce che un migliore stato di salute generale, misurato tramite l’esame, è correlato a una fertilità più bassa.
• Educazione (Education): Un aumento nel livello di educazione è correlato a una riduzione della fertilità. Le regioni con un livello di istruzione più elevato tendono ad avere una fertilità più bassa.

[Torna all’ Indice]

Interpretazione Geometrica

L’interpretazione geometrica della regressione lineare è un approccio concettuale che utilizza uno spazio tridimensionale (o superiore) per rappresentare visivamente il modello di regressione. In questo spazio, ogni punto rappresenta un’osservazione nel dataset, e un piano (o iperpiano) rappresenta il modello di regressione. L’obiettivo è trovare il piano (o iperpiano) che minimizza la somma dei quadrati delle distanze verticali tra i punti dati e il piano (o iperpiano). Questo fornisce una visualizzazione intuitiva di come i coefficienti del modello vengono stimati per ottenere la migliore “ajustement” ai dati, minimizzando gli errori residui. L’interpretazione geometrica aiuta a comprendere i principi fondamentali della regressione lineare e può essere applicata a problemi più complessi con più variabili indipendenti.

library(knitr)
library(plotly)

Loading required package: ggplot2
Registered S3 method overwritten by 'data.table':
  method           from
  print.data.table     
Registered S3 method overwritten by 'htmlwidgets':
  method           from         
  print.htmlwidget tools:rstudio

Attaching package: ‘plotly’

The following object is masked from ‘package:ggplot2’:

    last_plot

The following object is masked from ‘package:stats’:

    filter

The following object is masked from ‘package:graphics’:

    layout

# Genera dati casuali
set.seed(123)
n <- 50
X1 <- rnorm(n)
X2 <- rnorm(n)
Y <- 2 * X1 + 3 * X2 + rnorm(n)

# Crea un dataframe con le variabili
data <- data.frame(X1, X2, Y)

# Adatta il modello di regressione lineare
model <- lm(Y ~ X1 + X2, data = data)
summary(model)

Call:
lm(formula = Y ~ X1 + X2, data = data)

Residuals:
     Min       1Q   Median       3Q      Max 
-1.88137 -0.74056 -0.06374  0.52516  2.27045 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  -0.2299     0.1431  -1.607    0.115    
X1            2.0252     0.1540  13.150   <2e-16 ***
X2            2.8304     0.1575  17.973   <2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 0.9975 on 47 degrees of freedom
Multiple R-squared:  0.9107,    Adjusted R-squared:  0.9069 
F-statistic: 239.8 on 2 and 47 DF,  p-value: < 2.2e-16

# Creazione di una griglia di punti
x1_range <- seq(min(X1), max(X1), length = 20)
x2_range <- seq(min(X2), max(X2), length = 20)
grid <- expand.grid(X1 = x1_range, X2 = x2_range)

# Calcolo delle previsioni del modello sulla griglia
grid$Y_pred <- predict(model, newdata = grid)

# Creazione del plot 3D con plotly
 plot_ly(data, x = ~X1, y = ~X2, z = ~Y, type = "scatter3d", mode = "markers", marker = list(size = 5, color = "blue")) %>%
   add_surface(
     x = x1_range,
     y = x2_range,
     z = matrix(grid$Y_pred, nrow = length(x1_range), ncol = length(x2_range), byrow = TRUE),
     colors = "red",
     opacity = 0.7
   ) %>%
   layout(scene = list(xaxis = list(title = "X1"), yaxis = list(title = "X2"), zaxis = list(title = "Y")))

NA

L’interpretazione geometrica ci consente di vedere come il piano di regressione si adatta ai dati nello spazio tridimensionale e come i coefficienti stimati influenzano la posizione e l’inclinazione del piano rispetto ai dati osservati. Questo fornisce una visualizzazione intuitiva della relazione tra le variabili indipendenti e dipendenti nel contesto della regressione lineare.

[Torna all’ Indice]

Distribuzione F & Anova Table

L’ANOVA valuta globalmente se almeno una delle variabili indipendenti ha un effetto significativo sulla variabile dipendente, fornisce una statistica F e il relativo p-value. Un p-value basso suggerisce che almeno una delle variabili indipendenti è significativa nel modello. L’ANOVA fornisce quindi una visione complessiva della significatività del modello nel suo complesso.

L’analisi della varianza (ANOVA) e i “signif. codes” nel summary del modello forniscono informazioni simili, ma si concentrano su aspetti diversi dell’analisi.

I “signif. codes” nel summary del modello forniscono una valutazione variabile per variabile, indicando la significatività statistica di ciascun coefficiente. Utilizza asterischi (*) o altri simboli per indicare il livello di significatività, ad esempio, “***” potrebbe indicare un livello di significatività molto elevato (p-value molto basso), mentre ” ” (spazio) potrebbe indicare non significativo. Questa parte del summary fornisce una visione più dettagliata sulla significatività di ciascuna variabile indipendente separatamente.

L’ANOVA valuta la significatività del modello nel suo insieme, mentre i “signif. codes” nel summary forniscono una visione dettagliata della significatività di ciascuna variabile indipendente.

Esempio 1:

# Creiamo un dataset fittizio
set.seed(123)
data <- data.frame(
  Gruppo = rep(c("A", "B", "C"), each = 20),
  Punteggio = rnorm(60, mean = c(70, 75, 80), sd = 5)
)

# Eseguiamo l'ANOVA
anova_result <- aov(Punteggio ~ Gruppo, data = data)

# Visualizziamo la tabella ANOVA
summary(anova_result)

            Df Sum Sq Mean Sq F value Pr(>F)
Gruppo       2   11.3    5.65   0.143  0.867
Residuals   57 2255.9   39.58               

In questo esempio, eseguiamo un’ANOVA a un fattore per valutare le differenze nei punteggi tra i gruppi A, B e C. Il rapporto F e il valore p ci permettono di determinare se le differenze tra i gruppi sono statisticamente significative.

I risultati suggeriscono che non ci sono differenze statisticamente significative tra i gruppi, ossia il variare del gruppo non influenza significativamente la variabile dipendente. La varianza tra i gruppi è molto piccola rispetto alla varianza all’interno dei gruppi, e il test F non è significativo. Questo può indicare che i gruppi sono simili tra loro per quanto riguarda la variabile in studio.

Esempio 2:

# Creiamo un dataset fittizio
set.seed(123)
data <- data.frame(
  Genere = rep(c("Maschio", "Femmina"), each = 50),
  Trattamento = rep(c("A", "B"), times = 50),
  Punteggio = rnorm(100, mean = c(75, 80), sd = 5)
)

# Eseguiamo l'ANOVA a due fattori
anova_result <- aov(Punteggio ~ Genere * Trattamento, data = data)

# Visualizziamo la tabella ANOVA
summary(anova_result)

                   Df Sum Sq Mean Sq F value   Pr(>F)    
Genere              1    7.8     7.8   0.367    0.546    
Trattamento         1  575.8   575.8  26.953 1.16e-06 ***
Genere:Trattamento  1    2.6     2.6   0.121    0.729    
Residuals          96 2050.8    21.4                     
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

In questo esempio, eseguiamo un’ANOVA a due fattori per esaminare le differenze nei punteggi in base al genere e al trattamento.

Dalla tabella dell’ANOVA fornita, possiamo trarre le seguenti conclusioni:

• Il fattore “Trattamento” ha un forte effetto sulla variabile dipendente, con un valore di F elevato e un p-value molto basso, indicando che le differenze tra i trattamenti sono statisticamente significative (***).
• Il fattore “Genere” non ha un effetto significativo sulla variabile dipendente, con un valore di F basso e un p-value elevato.
• L’interazione tra “Genere” e “Trattamento” non ha un effetto significativo sulla variabile dipendente, con un valore di F e un p-value non significativi.

I risultati dell’ANOVA suggeriscono che il “Trattamento” è il principale driver delle differenze osservate nella variabile dipendente, mentre il “Genere” e l’interazione tra “Genere” e “Trattamento” non sembrano avere un effetto significativo.

• Un valore F maggiore di 1 suggerisce che i parametri o i fattori sono significativi, poiché la varianza spiegata è maggiore della varianza non spiegata.
• Un valore F vicino a 1 indica che il modello non spiega in modo significativo la variabilità nei dati.
• Il valore p associato all’F-value fornisce la probabilità che i risultati osservati siano dovuti al caso. Un valore p basso (di solito inferiore a 0.05) indica una significatività elevata, mentre un valore p alto suggerisce una mancanza di significatività.

Sistema di ipotesi:

Nell’analisi statistica in cui si calcola un valore F, ci sono due ipotesi principali: l’ipotesi nulla (H0) e l’ipotesi alternativa (H1).

Ipotesi Nulla (H0): L’ipotesi nulla afferma che non ci sono differenze significative tra i gruppi o i fattori considerati. In altre parole, l’ipotesi nulla sostiene che i parametri del modello o i fattori non hanno un effetto significativo sul risultato o che le differenze osservate sono casuali.

Ipotesi Alternativa (H1 o HA): L’ipotesi alternativa è il contrario dell’ipotesi nulla. Sostiene che ci sono differenze significative tra i gruppi o i fattori considerati, e che le differenze osservate non sono casuali, ma sono dovute a un effetto significativo dei parametri del modello o dei fattori.

Continuiamo l’esempio di prima:

model <- lm(Punteggio ~ Genere * Trattamento, data = data)
summary(model)

Call:
lm(formula = Punteggio ~ Genere * Trattamento, data = data)

Residuals:
     Min       1Q   Median       3Q      Max 
-12.3383  -3.1482  -0.1363   3.1717  10.9341 

Coefficients:
                           Estimate Std. Error t value Pr(>|t|)    
(Intercept)                 75.6716     0.9244  81.861  < 2e-16 ***
GenereMaschio               -0.2383     1.3073  -0.182 0.855741    
TrattamentoB                 5.1208     1.3073   3.917 0.000168 ***
GenereMaschio:TrattamentoB  -0.6434     1.8488  -0.348 0.728578    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 4.622 on 96 degrees of freedom
Multiple R-squared:  0.2223,    Adjusted R-squared:  0.198 
F-statistic: 9.147 on 3 and 96 DF,  p-value: 2.206e-05

# Estraiamo il valore p
p_value <- anova_summary[[1]][["Pr(>F)"]][3]  # Usiamo [3] per estrarre il valore relativo all'interazione


# Scegliamo un livello di significatività (alpha)
alpha <- 0.05

# Valutiamo se rifiutare l'ipotesi nulla
if (p_value < alpha) {
  cat("Rifiutiamo l'ipotesi nulla. Ci sono differenze significative tra i gruppi.\n")
} else {
  cat("Non rifiutiamo l'ipotesi nulla. Non ci sono differenze significative tra i gruppi.\n")
}

Non rifiutiamo l'ipotesi nulla. Non ci sono differenze significative tra i gruppi.

La F-statistic nel summary di un modello lineare (lm) rappresenta la statistica del test F per l’intero modello. Questo test verifica se c’è almeno una variabile indipendente nel modello che è significativamente associata alla variabile dipendente. In altre parole, valuta l’ipotesi nulla che tutti i coefficienti delle variabili indipendenti nel modello siano uguali a zero (cioè che non ci siano effetti).

In breve la F-statistic suggerisce che almeno una delle interazioni tra “Genere” e “Trattamento” o almeno una delle principali effetti è significativa nel modello.

[Torna all’ Indice]

Modelli Nidificati

Nei modelli statistici, un “nested model” si verifica quando un modello più complesso o generale può essere suddiviso o semplificato in un modello più semplice o specifico. Il modello più semplice è considerato “nidificato” all’interno del modello più complesso, poiché contiene un sottoinsieme di parametri o vincoli del modello più generale.

Nel contesto della regressione, i modelli nidificati sono spesso utilizzati per testare l’aggiunta di variabili indipendenti al modello al fine di valutare se le variabili aggiuntive migliorano significativamente la capacità di previsione o spiegazione del modello. I modelli nidificati sono anche utilizzati in contesti come l’analisi della varianza (ANOVA), l’analisi della devianza nei modelli lineari generalizzati (che vedremo più avanti) e altre procedure statistiche.

Esempio di Modelli di Regressione Nidificati:

Supponiamo di voler creare un modello di regressione per prevedere il reddito di una persona basandoci su quattro variabili indipendenti: età, istruzione, esperienza lavorativa e genere. Il modello completo potrebbe essere:

Modello Completo

\[ Reddito = \beta_0 + \beta_1 \cdot Età + \beta_2 \cdot Istruzione + \beta_3 \cdot Esperienza + \beta_4 \cdot Genere \] Tuttavia, potremmo essere interessati a valutare se l’aggiunta della variabile “genere” migliora significativamente la capacità predittiva del modello. In tal caso, il modello senza “genere” è nidificato all’interno del modello completo:

Modello Nidificato

\[ Reddito' = \beta_0 + \beta_1 \cdot Età + \beta_2 \cdot Istruzione + \beta_3 \cdot Esperienza \]

In questo esempio, il Modello 1 è il modello completo e il Modello 2 è il modello nidificato senza il parametro per “genere”. Valutiamo tramite un ANOVA quale dei due modelli sia considerato più significativo.

# Creiamo dati fittizi
set.seed(123)
n <- 100
eta <- rnorm(n, mean = 35, sd = 5)
istruzione <- rnorm(n, mean = 12, sd = 2)
esperienza <- rnorm(n, mean = 10, sd = 3)
genere <- sample(c("Maschio", "Femmina"), n, replace = TRUE)
reddito <- 20 + 2 * eta + 3 * istruzione + 5 * esperienza + ifelse(genere == "Maschio", 4, 0) + rnorm(n, mean = 0, sd = 5)

# Creiamo un dataframe con i dati
data <- data.frame(eta, istruzione, esperienza, genere, reddito)

# Modello completo
modello_completo <- lm(reddito ~ eta + istruzione + esperienza + genere, data = data)

# Modello nidificato senza "genere"
modello_nidificato <- lm(reddito ~ eta + istruzione + esperienza, data = data)

# Test F per confrontare i modelli
anova_result <- anova(modello_nidificato, modello_completo)

# Visualizziamo la tabella ANOVA
print(anova_result)

Analysis of Variance Table

Model 1: reddito ~ eta + istruzione + esperienza
Model 2: reddito ~ eta + istruzione + esperienza + genere
  Res.Df    RSS Df Sum of Sq      F    Pr(>F)    
1     96 3137.3                                  
2     95 2606.6  1    530.69 19.342 2.852e-05 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

# Scegliamo un livello di significatività (alpha)
alpha <- 0.05

# Valutiamo se rifiutare l'ipotesi nulla
if (anova_result[2, "Pr(>F)"] < alpha) {
  cat("La rimozione di 'genere' non migliora significativamente il modello.\n")
} else {
  cat("La rimozione di 'genere' migliora significativamente il modello.\n")
}

La rimozione di 'genere' non migliora significativamente il modello.

# R Squared dei due modelli
print("Modello completo ") 

[1] "Modello completo "

summary(modello_completo)$r.squared

[1] 0.9155659

print("Modello annidato ")

[1] "Modello annidato "

summary(modello_nidificato)$r.squared

[1] 0.8983755

[Torna all’ Indice]

Selezione delle variabili

La “variable selection” è un processo attraverso il quale si scelgono le variabili più rilevanti da includere in un modello statistico. Questo processo è utile per semplificare i modelli, migliorare la capacità predittiva e la comprensione dei dati, ridurre l’overfitting e aumentare l’efficienza computazionale.

Un metodo comune per la selezione delle variabili in R coinvolge l’utilizzo dell’Information Criterion (Criterio d’Informazione) di Akaike (AIC) insieme alla funzione step().

AIC (Akaike’s Information Criterion):

Il Criterio d’Informazione di Akaike (AIC) è una metrica che misura la qualità di un modello statistico. L’obiettivo dell’AIC è trovare il miglior compromesso tra la bontà di adattamento del modello ai dati e la sua complessità. L’AIC tiene conto della funzione di verosimiglianza del modello e penalizza i modelli con un numero maggiore di parametri. L’AIC è definito come:

\[ AIC = -2logLikelihood + 2k \]

Dove:

• “log-likelihood” è il logaritmo della funzione di verosimiglianza del modello.
• “k” è il numero di parametri stimati nel modello. Un valore AIC più basso indica un modello migliore, in quanto indica un migliore adattamento ai dati con meno complessità.

Funzione step():

La funzione step() in R è utilizzata per effettuare la selezione delle variabili basata su criteri come l’AIC. Consente di confrontare e selezionare i modelli in modo automatico aggiungendo o rimuovendo variabili dal modello, fino a trovare il modello con l’AIC più basso. La sintassi di base della funzione step() è la seguente:

#step(modello_iniziale, direction = "both", scope = list(lower = modello_minimo, upper = modello_massimo))

• modello_iniziale è il modello di partenza che desideri semplificare o migliorare.
• direction può essere “forward”, “backward”, o “both” e specifica se aggiungere, rimuovere o entrambi i tipi di variabili durante la selezione.
• scope specifica l’intervallo dei modelli da considerare durante la selezione. Il - “modello_minimo” rappresenta il modello più semplice possibile (ad esempio, un modello con solo l’intercetta), mentre il “modello_massimo” rappresenta il modello più complesso (il modello completo con tutte le variabili).

# Carica il dataset di esempio
data(mtcars)

# Crea un modello lineare iniziale
all <- lm(mpg ~ ., data = mtcars)

# Esegui la selezione delle variabili basata su AIC
best <- step(all, direction = "backward")

Start:  AIC=70.9
mpg ~ cyl + disp + hp + drat + wt + qsec + vs + am + gear + carb

       Df Sum of Sq    RSS    AIC
- cyl   1    0.0799 147.57 68.915
- vs    1    0.1601 147.66 68.932
- carb  1    0.4067 147.90 68.986
- gear  1    1.3531 148.85 69.190
- drat  1    1.6270 149.12 69.249
- disp  1    3.9167 151.41 69.736
- hp    1    6.8399 154.33 70.348
- qsec  1    8.8641 156.36 70.765
<none>              147.49 70.898
- am    1   10.5467 158.04 71.108
- wt    1   27.0144 174.51 74.280

Step:  AIC=68.92
mpg ~ disp + hp + drat + wt + qsec + vs + am + gear + carb

       Df Sum of Sq    RSS    AIC
- vs    1    0.2685 147.84 66.973
- carb  1    0.5201 148.09 67.028
- gear  1    1.8211 149.40 67.308
- drat  1    1.9826 149.56 67.342
- disp  1    3.9009 151.47 67.750
- hp    1    7.3632 154.94 68.473
<none>              147.57 68.915
- qsec  1   10.0933 157.67 69.032
- am    1   11.8359 159.41 69.384
- wt    1   27.0280 174.60 72.297

Step:  AIC=66.97
mpg ~ disp + hp + drat + wt + qsec + am + gear + carb

       Df Sum of Sq    RSS    AIC
- carb  1    0.6855 148.53 65.121
- gear  1    2.1437 149.99 65.434
- drat  1    2.2139 150.06 65.449
- disp  1    3.6467 151.49 65.753
- hp    1    7.1060 154.95 66.475
<none>              147.84 66.973
- am    1   11.5694 159.41 67.384
- qsec  1   15.6830 163.53 68.200
- wt    1   27.3799 175.22 70.410

Step:  AIC=65.12
mpg ~ disp + hp + drat + wt + qsec + am + gear

       Df Sum of Sq    RSS    AIC
- gear  1     1.565 150.09 63.457
- drat  1     1.932 150.46 63.535
<none>              148.53 65.121
- disp  1    10.110 158.64 65.229
- am    1    12.323 160.85 65.672
- hp    1    14.826 163.35 66.166
- qsec  1    26.408 174.94 68.358
- wt    1    69.127 217.66 75.350

Step:  AIC=63.46
mpg ~ disp + hp + drat + wt + qsec + am

       Df Sum of Sq    RSS    AIC
- drat  1     3.345 153.44 62.162
- disp  1     8.545 158.64 63.229
<none>              150.09 63.457
- hp    1    13.285 163.38 64.171
- am    1    20.036 170.13 65.466
- qsec  1    25.574 175.67 66.491
- wt    1    67.572 217.66 73.351

Step:  AIC=62.16
mpg ~ disp + hp + wt + qsec + am

       Df Sum of Sq    RSS    AIC
- disp  1     6.629 160.07 61.515
<none>              153.44 62.162
- hp    1    12.572 166.01 62.682
- qsec  1    26.470 179.91 65.255
- am    1    32.198 185.63 66.258
- wt    1    69.043 222.48 72.051

Step:  AIC=61.52
mpg ~ hp + wt + qsec + am

       Df Sum of Sq    RSS    AIC
- hp    1     9.219 169.29 61.307
<none>              160.07 61.515
- qsec  1    20.225 180.29 63.323
- am    1    25.993 186.06 64.331
- wt    1    78.494 238.56 72.284

Step:  AIC=61.31
mpg ~ wt + qsec + am

       Df Sum of Sq    RSS    AIC
<none>              169.29 61.307
- am    1    26.178 195.46 63.908
- qsec  1   109.034 278.32 75.217
- wt    1   183.347 352.63 82.790

In questo esempio, partiamo da un modello lineare completo che utilizza tutte le variabili di mtcars, e poi utilizziamo step() per eseguire la selezione delle variabili basata su AIC. Alla fine, otteniamo il modello con l’AIC più basso, che dovrebbe essere una versione semplificata del modello iniziale con solo le variabili più rilevanti.

La “variable selection” utilizzando AIC e step() è un potente strumento per migliorare la qualità e l’interpretabilità dei modelli statistici, in particolare quando si hanno molti potenziali predittori.

[Torna all’ Indice]

Predizioni Categoriche

Nell’analisi statistica, i “categorical predictors” sono variabili che rappresentano categorie o gruppi distinti anziché valori numerici. Queste variabili sono anche conosciute come variabili qualitative o fattori. Ad esempio, il genere (maschio/femmina), il livello di istruzione (scuola elementare, scuola media, laurea), o il tipo di prodotto (A, B, C) sono esempi di predittori categorici. Quando si utilizzano predittori categorici in un modello statistico, è importante considerare come gestire e interpretare questi dati.

Una considerazione fondamentale è come rappresentare le variabili categoriche nel modello. Solitamente, vengono utilizzate delle variabili dummy (variabili indicatrici) per rappresentare le categorie. Ad esempio, nel caso del genere (maschio/femmina), potrebbero essere create due variabili dummy, una per il maschio e una per la femmina. Queste variabili dummy prendono il valore 1 o 0 a seconda dell’appartenenza alla categoria. Questo approccio consente al modello di catturare l’effetto della categoria sulla variabile dipendente.

Oltre alla rappresentazione delle variabili categoriche, è importante considerare le interazioni tra i predittori. Le interazioni si verificano quando l’effetto di una variabile categorica sul risultato dipende da un’altra variabile. Ad esempio, l’effetto del livello di istruzione sul reddito potrebbe variare in base al genere. In questo caso, c’è un’interazione tra il livello di istruzione e il genere.

Per esaminare le interazioni tra predittori categorici, è possibile utilizzare l’analisi della varianza (ANOVA) o i modelli lineari generalizzati (che vedremo più avanti). Le interazioni possono fornire informazioni preziose sull’influenza combinata delle variabili categoriche sul risultato.

# Creiamo dati fittizi
set.seed(123)
n <- 100
genere <- sample(c("Maschio", "Femmina"), n, replace = TRUE)
istruzione <- rep(c("Elementare", "Media", "Laurea"), length.out = n )
reddito <- 30 + ifelse(genere == "Maschio", 5, 0) + ifelse(istruzione == "Laurea", 10, 0) + rnorm(n, mean = 0, sd = 5)

# Creiamo un dataframe con i dati
data <- data.frame(genere, istruzione, reddito)

# Modello con interazione tra genere e istruzione
modello <- lm(reddito ~ genere * istruzione, data = data)

# Visualizziamo i risultati
summary(modello)

Call:
lm(formula = reddito ~ genere * istruzione, data = data)

Residuals:
    Min      1Q  Median      3Q     Max 
-9.3888 -3.0121 -0.6269  2.6272 11.0710 

Coefficients:
                               Estimate Std. Error t value Pr(>|t|)    
(Intercept)                     27.8429     1.2794  21.763  < 2e-16 ***
genereMaschio                    5.8726     1.6681   3.521 0.000666 ***
istruzioneLaurea                12.0839     1.8438   6.554 2.98e-09 ***
istruzioneMedia                  3.3778     1.7519   1.928 0.056860 .  
genereMaschio:istruzioneLaurea  -1.3699     2.3856  -0.574 0.567186    
genereMaschio:istruzioneMedia   -0.8084     2.3586  -0.343 0.732535    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 4.787 on 94 degrees of freedom
Multiple R-squared:  0.5839,    Adjusted R-squared:  0.5618 
F-statistic: 26.38 on 5 and 94 DF,  p-value: < 2.2e-16

anova(modello, test = "chi")

Analysis of Variance Table

Response: reddito
                  Df  Sum Sq Mean Sq F value    Pr(>F)    
genere             1  738.22  738.22 32.2149 1.525e-07 ***
istruzione         2 2277.06 1138.53 49.6837 1.892e-15 ***
genere:istruzione  2    7.66    3.83  0.1672    0.8463    
Residuals         94 2154.06   22.92                      
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

# Grafico reddito per soli maschi e sole femmine
ggplot(data = data, aes(x = istruzione, y = reddito, fill = genere)) +
  geom_boxplot() +
  labs(x = "Istruzione", y = "Reddito")

In questo esempio, stiamo creando dati fittizi con due predittori categorici: “genere” e “istruzione”. Il modello lineare include un’interazione tra questi due predittori. La tabella dei risultati summary(modello) mostra come i predittori categorici e l’interazione influenzano il reddito.

Possiamo concludere che il genere e il livello di istruzione hanno un effetto significativo sul reddito, mentre le interazioni tra genere e istruzione non sono significative in questo modello. Il modello nel suo complesso è significativo e in grado di spiegare una parte della variazione nel reddito.

Fattori con più di due categorie:

Quando si affrontano fattori con più di due livelli (categorie), è necessario considerare come gestire queste variabili nel modello. In generale, un fattore con k livelli richiede la creazione di k-1 variabili dummy per evitare la “dummy variable trap”. Questo si verifica quando le variabili dummy sono linearmente dipendenti e possono portare a problemi di multicollinearità.

Ad esempio, se abbiamo una variabile “colore” con tre livelli (rosso, verde, blu), dovremmo creare due variabili dummy per rappresentarla. Una rappresenterà il rosso e l’altra il verde. Se entrambe le variabili dummy sono uguali a 0, ciò significa che il colore è blu. Questo evita la trap della variabile dummy.

# Creiamo dati fittizi
set.seed(123)
n <- 100
colore <- rep(c("Rosso", "Verde", "Blu"), length.out = n )
voto <- rnorm(n, mean = 50, sd = 10)

# Creiamo un dataframe con i dati
data <- data.frame(colore, voto)

# Modello con un fattore con più di due livelli
modello <- lm(voto ~ colore, data = data)

# Visualizziamo i risultati
summary(modello)

Call:
lm(formula = voto ~ colore, data = data)

Residuals:
     Min       1Q   Median       3Q      Max 
-23.6349  -6.1247  -0.4198   5.7870  21.9547 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  50.5432     1.5932  31.725   <2e-16 ***
coloreRosso   1.8459     2.2364   0.825    0.411    
coloreVerde  -0.8084     2.2531  -0.359    0.721    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 9.152 on 97 degrees of freedom
Multiple R-squared:  0.01508,   Adjusted R-squared:  -0.005228 
F-statistic: 0.7426 on 2 and 97 DF,  p-value: 0.4786

I risultati indicano che il colore del prodotto (rosso o verde) non ha un impatto significativo sul voto. L’intercetta, che rappresenta il colore “Blu,” è significativa, ma il modello nel suo insieme non è molto efficace nel spiegare la variazione nei voti.

[Torna all’ Indice]

Controllo del Modello

Il Model Checking è una fase cruciale nell’analisi statistica, specialmente quando si adotta un modello di regressione. Durante questa fase, si valuta se il modello soddisfa le principali assunzioni dei modelli lineari. Le quattro assunzioni principali da verificare sono:

• Linearity (Linearità): Questa assunzione afferma che la risposta (variabile dipendente) può essere scritta come una combinazione lineare delle variabili predittive (variabili indipendenti). In altre parole, il modello dovrebbe essere in grado di catturare il rapporto tra le variabili in modo lineare, con un certo grado di rumore residuo. La linearità può essere verificata attraverso grafici di dispersione o grafici residui.
• Independence (Indipendenza): Questa assunzione richiede che gli errori (residui) del modello siano indipendenti l’uno dall’altro. Ciò significa che il valore di errore per un’osservazione non è influenzato dal valore di errore per un’altra osservazione. L’indipendenza può essere verificata osservando i grafici dei residui in sequenza temporale o spaziale, a seconda del contesto.
• Normality (Normalità): L’assunzione di normalità richiede che i residui del modello seguano una distribuzione normale. Questo è importante perché molte procedure statistiche si basano sull’ipotesi di normalità dei residui. La normalità può essere verificata tramite grafici quantile-quantile (QQ plot) o istogrammi dei residui.
• Equal Variance (Varianza Uniforme): Questa assunzione, chiamata anche omoschedasticità, richiede che la varianza dei residui sia costante in tutti i livelli delle variabili predittive. In altre parole, non dovrebbe esserci alcun modello discernibile nella varianza dei residui. La varianza uniforme può essere verificata osservando i grafici dei residui rispetto ai valori predetti.

Per verificare queste assunzioni, i Residuals-based displays (grafici basati sui residui) sono spesso utilizzati. Questi includono:

• Scatterplot dei residui: Un grafico dei residui contro i valori previsti o le variabili predittive. Questo può rivelare se c’è una struttura non lineare nei residui.
• Grafico di sequenza temporale dei residui: Utilizzato quando i dati sono raccolti nel tempo, questo grafico può rivelare dipendenze temporali nei residui.
• QQ-plot (Quantile-Quantile plot): Questo grafico confronta i quantili dei residui con quelli di una distribuzione normale. Se i punti del grafico seguono una linea retta, i residui sono approssimativamente normali.
• Istogramma dei residui: Un istogramma dei residui può dare un’idea della loro distribuzione e normalità.

Rispettare queste assunzioni è importante per garantire che le stime del modello siano affidabili e che le conclusioni siano valide. Se una o più di queste assunzioni non sono soddisfatte, potrebbero essere necessarie correzioni al modello o ai dati stessi.

[Torna all’ Indice]

Transformazioni

Le trasformazioni sono una tecnica utilizzata nella modellazione statistica per modificare le relazioni tra variabili al fine di soddisfare meglio le assunzioni del modello. Le trasformazioni possono essere utili quando le relazioni tra le variabili non sono lineari o quando le assunzioni di omoschedasticità o normalità dei residui non sono soddisfatte. Di seguito, affrontiamo i seguenti argomenti relativi alle trasformazioni:

• Variance Stabilizing Transformations (Trasformazioni per Stabilizzare la Varianza): In alcuni casi, la varianza dei dati può variare in modo non costante con il cambiare del valore medio. Questo fenomeno è noto come eteroschedasticità. Le trasformazioni possono essere utilizzate per stabilizzare la varianza, rendendo la relazione tra il valore medio e la varianza più costante. Un esempio comune è la trasformazione di Box-Cox.

• Box-Cox Transform: La trasformazione di Box-Cox è una tecnica utilizzata per stabilizzare la varianza e rendere i dati approssimativamente normali. È definita come:

\[ y(\lambda) = \begin{cases} \frac{(y^\lambda - 1)}{\lambda} & \text{se } \lambda \neq 0 \\ \log(y) & \text{se } \lambda = 0 \end{cases} \]

Dove y sono i dati originali e λ è il parametro di trasformazione. È possibile calcolare il valore ottimale di λ che massimizza la normalità dei dati.

library(ggplot2)
library(MASS)
library(gridExtra)

# Genera dati casuali
set.seed(123)
data <- data.frame(y = rgamma(100, shape = 2, scale = 1))

# Applica la trasformazione di Box-Cox
result <- boxcox(y ~ 1, data = data)

lambda <- result$x[which.max(result$y)]
transformed_data <- if (lambda == 0) log(data$y) else ((data$y^lambda - 1) / lambda)

# Visualizza il valore ottimale di lambda
cat("Valore ottimale di lambda: ", lambda, "\n")

Valore ottimale di lambda:  0.3838384 

# Crea un dataframe con i dati originali e trasformati
plot_data <- data.frame(Original = data$y, Transformed = transformed_data)

# Plotta i dati originali
plot_original <- ggplot(plot_data, aes(x = Original)) +
  geom_histogram(binwidth = 0.5, fill = "blue", alpha = 0.7) +
  labs(title = "Distribuzione dei dati originali")

# Plotta i dati trasformati
plot_transformed <- ggplot(plot_data, aes(x = Transformed)) +
  geom_histogram(binwidth = 0.1, fill = "green", alpha = 0.7) +
  labs(title = "Distribuzione dei dati trasformati")

# Mostra i grafici sulla stessa riga
grid.arrange(plot_original, plot_transformed, ncol = 2)

• Polynomials (Polinomi): Le trasformazioni polinomiali consentono di modellare relazioni non lineari tra variabili. È possibile aggiungere termini polinomiali al modello di regressione per catturare curve o relazioni più complesse. Ad esempio, si possono utilizzare polinomi di secondo grado per modellare una relazione quadratica tra una variabile indipendente e la variabile dipendente. L’aggiunta di termini polinomiali può migliorare l’adattamento del modello ai dati, ma è importante evitare di aggiungere troppi termini polinomiali per evitare l’overfitting.

# Modello lineare con un termine polinomiale di secondo grado
model <- lm(y ~ x + I(x^2), data = data)

• Transformations of Predictor Variables (Trasformazioni delle Variabili Predittive): Le trasformazioni delle variabili predittive sono utilizzate per adattare i dati in modo che soddisfino meglio le assunzioni del modello. Queste trasformazioni coinvolgono la modifica delle variabili indipendenti piuttosto che della variabile dipendente. Possono essere utilizzate per rendere le relazioni tra le variabili più lineari o per stabilizzare la varianza. Ad esempio, è possibile applicare una trasformazione logaritmica o una radice quadrata a una variabile predittiva per renderla più lineare nei confronti della variabile dipendente.

# Creiamo dati fittizi
set.seed(123)
X <- rnorm(100, mean = 10, sd = 2)
Y <- 2 * X + rnorm(100, mean = 0, sd = 1)

# Creiamo un dataframe con i dati
data <- data.frame(X, Y)

# Modello lineare senza trasformazione
model_no_transform <- lm(Y ~ X, data = data)

# Visualizziamo il summary del modello senza trasformazione
summary(model_no_transform)

Call:
lm(formula = Y ~ X, data = data)

Residuals:
    Min      1Q  Median      3Q     Max 
-1.9073 -0.6835 -0.0875  0.5806  3.2904 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  0.15956    0.55265   0.289    0.773    
X            1.97376    0.05344  36.935   <2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 0.9707 on 98 degrees of freedom
Multiple R-squared:  0.933, Adjusted R-squared:  0.9323 
F-statistic:  1364 on 1 and 98 DF,  p-value: < 2.2e-16

# Trasformiamo la variabile X applicando il logaritmo
data$X_transformed <- log(data$X)

# Modello lineare con la variabile X trasformata
model_with_transform <- lm(Y ~ X_transformed, data = data)

# Visualizziamo il summary del modello con la variabile X trasformata
summary(model_with_transform)

Call:
lm(formula = Y ~ X_transformed, data = data)

Residuals:
    Min      1Q  Median      3Q     Max 
-1.9290 -0.7441 -0.1202  0.5388  3.3305 

Coefficients:
              Estimate Std. Error t value Pr(>|t|)    
(Intercept)   -24.1696     1.2961  -18.65   <2e-16 ***
X_transformed  19.2819     0.5608   34.38   <2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 1.037 on 98 degrees of freedom
Multiple R-squared:  0.9235,    Adjusted R-squared:  0.9227 
F-statistic:  1182 on 1 and 98 DF,  p-value: < 2.2e-16

Le trasformazioni sono strumenti potenti per adattare i modelli ai dati in modo più accurato quando le relazioni tra variabili non sono lineari o quando le assunzioni del modello non sono soddisfatte. Tuttavia, è importante scegliere con attenzione le trasformazioni per evitare il sovradattamento e garantire che i risultati siano interpretabili.

[Torna all’ Indice]

Multicollinearità

La multicollinearità si verifica quando due o più variabili indipendenti in un modello di regressione sono fortemente correlate tra loro. Questa correlazione tra le variabili indipendenti può rendere difficile l’interpretazione del modello e portare a stime poco affidabili dei coefficienti di regressione. La presenza di multicollinearità può causare un aumento del Variance Inflation Factor (VIF), una misura comune utilizzata per valutare la multicollinearità tra le variabili indipendenti in un modello di regressione. Un alto VIF per una variabile indica che quella variabile è fortemente correlata con le altre variabili indipendenti nel modello.

Variance Inflation Factor:

Il VIF di ciascuna variabile indipendente è calcolato come il rapporto della varianza dell’errore standard del coefficiente di regressione stimato per quella variabile rispetto alla varianza dell’errore standard se la variabile fosse stata completamente non correlata alle altre variabili indipendenti. In generale, un VIF superiore a 5 o 10 è spesso considerato un segno di multicollinearità significativa.

\[ VIF_i = (X^TX)_{i+1,i+1}^{-1}*ns^2_{X_i} \ \ oppure \ \ VIF_i = \frac{1}{1-R^2_i} \]

# Load the necessary library
library(car)

# Create a sample dataset with multiple predictor variables
set.seed(123)
data <- data.frame(
  X1 = rnorm(100),
  X2 = rnorm(100),
  X3 = rnorm(100),
  X4 = rnorm(100)
)

# Add a dependent variable (response)
data$Y <- 2 * data$X1 + 3 * data$X2 + 1.5 * data$X3 + rnorm(100)

# Fit a linear regression model
model <- lm(Y ~ X1 + X2 + X3 + X4, data = data)

# Calculate VIF
vif_values <- vif(model)

# Print the VIF values
vif_values

      X1       X2       X3       X4 
1.021515 1.004920 1.020348 1.006429 

Tutti i valori VIF sono vicini a 1, il che suggerisce che non c’è una forte multicollinearità tra le variabili predittive X1, X2, X3 e X4. Questo è un buon segno, poiché significa che le variabili non sono fortemente correlate tra loro.

Valori VIF più elevati indicano una multicollinearità più forte, e valori al di sopra di una certa soglia (ad esempio, VIF > 5) possono suggerire la necessità di affrontare la collinearità, ad esempio, rimuovendo una delle variabili predittive correlate.

[Torna all’ Indice]

Punti Influenti

I punti influenti si riferiscono a osservazioni nei dati che hanno un impatto significativo sui risultati di un’analisi statistica, come una regressione lineare. Questi punti possono influenzare la stima dei parametri del modello, i residui, i valori p, l’R-squared e altre statistiche di rilevanza. Ci sono diverse metriche utilizzate per identificare i punti influenti, tra cui Standardized Residuals, Studentized Residuals e Cook’s Distance.

• Standardized Residuals (Residui Standardizzati): Questi sono i residui divisi per la deviazione standard dei residui. Un residuo standardizzato è una misura di quanto un punto dato si discosti dalla linea di regressione in termini di deviazioni standard. I punti con residui standardizzati molto grandi (positivi o negativi) sono considerati influenti.
• Studentized Residuals (Residui Studentizzati): Questi sono i residui divisi per una stima della deviazione standard dell’errore residuo. I residui studentizzati sono utilizzati per valutare quanto un punto dato sia influente considerando l’effetto delle altre osservazioni nel dataset. I punti con residui studentizzati significativamente grandi in valore assoluto sono considerati influenti.
• Cook’s Distance (Distanza di Cook): Cook’s Distance è una metrica che combina l’effetto di un punto sui parametri del modello e il suo effetto sui residui. I punti con Cook’s Distance molto grandi sono considerati influenti. Cook’s Distance è spesso utilizzato per identificare punti che, se rimossi, avrebbero un impatto significativo sui risultati del modello.

Nel contesto della regressione, i punti influenti possono derivare da outlier nei dati, dati errati o punti che influenzano notevolmente la stima dei parametri. Identificare e trattare i punti influenti è importante per garantire che il modello di regressione sia affidabile e rappresenti accuratamente i dati. La rimozione di punti influenti può migliorare la bontà di adattamento del modello e l’accuratezza delle previsioni.

# Carichiamo il dataset di esempio
data(mtcars)

# Adattiamo un modello di regressione lineare
model <- lm(mpg ~ wt + hp, data = mtcars)

# Calcoliamo i residui standardizzati
standardized_residuals <- rstandard(model)

# Identifichiamo i punti influenti basati sui residui standardizzati
influential_points <- which(abs(standardized_residuals) > 2)

# Visualizziamo gli indici dei punti influenti
cat("Punti influenti basati sui residui standardizzati:", influential_points, "\n")

Punti influenti basati sui residui standardizzati: 17 18 20 

# Calcoliamo Cook's Distance
cook_distance <- cooks.distance(model)

# Identifichiamo i punti influenti basati su Cook's Distance
influential_points_cook <- which(cook_distance > 4 / length(cook_distance))

# Visualizziamo gli indici dei punti influenti basati su Cook's Distance
cat("Punti influenti basati su Cook's Distance:", influential_points_cook, "\n")

Punti influenti basati su Cook's Distance: 17 18 20 31 

par(mfrow = c(1,2))
# Grafico dei punti con evidenziazione dei punti influenti
plot(mtcars$wt, mtcars$mpg, pch = 16, main = "Scatter plot con Punti Influenti (Stan e Stud)",
     xlab = "Peso (wt)", ylab = "Miglia per gallone (mpg)")
points(mtcars$wt[influential_points], mtcars$mpg[influential_points], pch = 16, col = "red", cex = 1.5)

plot(mtcars$wt, mtcars$mpg, pch = 16, main = "Scatter plot con Punti Influenti (Cook)",
     xlab = "Peso (wt)", ylab = "Miglia per gallone (mpg)")
points(mtcars$wt[influential_points_cook], mtcars$mpg[influential_points_cook], pch = 16, col = "red", cex = 1.5)

Leverage:

Il “Leverage” è una misura utilizzata nell’analisi dei dati statistici per identificare punti influenti o osservazioni atipiche in un modello di regressione. Questa misura valuta quanto un’osservazione può influenzare i risultati del modello, in particolare i coefficienti di regressione. Il leverage è calcolato sulla base delle variabili predittive e può essere utilizzato per identificare le osservazioni che hanno un impatto significativo sul modello.

Identificazione dei punti influenti: I punti con un valore di leverage significativamente più alto degli altri sono quelli che possono influenzare notevolmente il modello. Puoi stabilire una soglia arbitraria o utilizzare metodi statistici per determinare quali punti sono influenti. Ad esempio, i punti con leverage superiore a 2 volte la media possono essere considerati influenti.

Esame dei punti influenti: Una volta identificati i punti influenti, è possibile esaminarli ulteriormente per determinare se sono effettivamente outliers o errori di misurazione. Potresti voler esaminare le osservazioni con elevate differenze tra i valori osservati e quelli previsti dal modello.

È importante notare che la rimozione dei punti influenti dovrebbe essere effettuata con cautela e solo se c’è una giustificazione valida. In alcuni casi, potresti scegliere di mantenere i punti influenti nel modello se ritieni che rappresentino informazioni significative o se hanno una spiegazione plausibile.

# Generiamo dati casuali
set.seed(123)
n <- 100
x <- rnorm(n)
y <- 2 * x + rnorm(n)

# Adattiamo un modello di regressione lineare
model <- lm(y ~ x)

# Calcoliamo i valori di leverage
leverage <- hatvalues(model)

# Identifichiamo i punti influenti
infl_points <- which(leverage > 2 * mean(leverage))

# Visualizziamo i punti influenti
print(infl_points)

 6 16 18 26 44 57 70 72 97 
 6 16 18 26 44 57 70 72 97 

par(mfrow = c(1,2))
# Plot del grafico
plot(x, y)
points(x[infl_points], y[infl_points], col = "blue", pch = 19)

# Plot del Leverage
plot(x, leverage)
points(x[infl_points], leverage[infl_points], col = "blue", pch = 19)

Nota che nell’esempio abbiamo utilizzato una soglia di leverage arbitraria (2 volte la media) per identificare i punti influenti. In un’applicazione pratica, è consigliabile considerare la soglia in base al contesto del problema e all’analisi dei dati.

[Torna all’ Indice]

Modelli Lineari Generalizzati

I GLM estendono il framework della regressione lineare per gestire un’ampia gamma di distribuzioni dei dati e tipologie di risposte. A differenza della regressione lineare tradizionale, i GLM possono accomodare distribuzioni di errori non normali e modellare relazioni tra predittori e risposte attraverso una funzione di collegamento.

Un GLM è caratterizzato da tre componenti principali:

1. Componente Casuale (Distribuzione): La variabile di risposta \(Y\) segue una distribuzione di probabilità dalla famiglia esponenziale, che include distribuzioni comuni come normale, binomiale e di Poisson.

2. Componente Sistematica (Predittore Lineare): La relazione tra i predittori e il valore atteso della risposta è espressa attraverso un predittore lineare (\(\eta\)). Il predittore lineare è una combinazione dei predittori, ognuno moltiplicato per un parametro, e si collega alla media della risposta attraverso una funzione di collegamento.

3. Funzione di Collegamento: La funzione di collegamento (\(g(\mu)\)) stabilisce il collegamento tra il predittore lineare e la media della risposta. Trasforma la scala della variabile di risposta e assicura che il predittore lineare copra l’intera linea reale. Le funzioni di collegamento comuni includono logit, probit e identità.

La forma generale di un GLM può essere rappresentata come segue:

\[ g(\mu) = X\beta \]

Dove:

• \(g(\mu)\) è la funzione di collegamento.
• \(\mu\) è il valore atteso della variabile di risposta.
• \(X\) è la matrice dei predittori.
• \(\beta\) è il vettore dei coefficienti.

Esempi:

1. Regressione Logistica Binaria:

   • Distribuzione: Binomiale
   • Funzione di Collegamento: Logit (log-rapporti di probabilità)
   • Equazione: \(\text{logit}(\mu) = X\beta\)

   Descrizione: La funzione di collegamento logit trasforma la probabilità di successo (\(\mu\)) in un predittore lineare. In questo caso, il modello logistic descrive come la log-odds della probabilità di successo sia lineare rispetto ai predittori.

2. Regressione di Poisson:

   • Distribuzione: Poisson
   • Funzione di Collegamento: Log
   • Equazione: \(\log(\mu) = X\beta\)

   Descrizione: Con la funzione di collegamento logaritmico, il modello di Poisson può gestire dati di conteggio, poiché connette il logaritmo naturale del valore atteso (\(\mu\)) a un predittore lineare.

3. Regressione Gamma:

   • Distribuzione: Gamma
   • Funzione di Collegamento: Inverso
   • Equazione: \(\frac{1}{\mu} = X\beta\)

   Descrizione: La funzione di collegamento inversa in un modello gamma è appropriata quando si modellano variabili con distribuzioni a coda pesante. Collega l’inverso del valore atteso (\(\mu\)) a un predittore lineare.

Devianza

La devianza è una misura della discrepanza tra il modello statistico e i dati osservati nei GLM. In generale, la devianza è utilizzata per confrontare modelli alternativi e valutare quanto bene un modello si adatta ai dati. Nel contesto dei GLM, la devianza è particolarmente significativa perché tiene conto delle specifiche distribuzioni delle variabili di risposta.

La devianza si calcola confrontando il modello fitted (previsto) con un modello null, spesso noto come modello null di saturazione. Il modello null rappresenta l’ipotesi che tutti i parametri del modello siano uguali a zero, indicando l’assenza di effetti predittori. La devianza è data dalla seguente formula:

\[ D = 2 \times \left( \ell(\hat{\beta}) - \ell(\beta_0) \right) \]

Dove:

• \(\ell(\hat{\beta})\) è il log-likelihood del modello fitted.
• \(\ell(\beta_0)\) è il log-likelihood del modello null.
• \(D\) è la devianza.

Poiché i GLM utilizzano la famiglia esponenziale di distribuzioni, la devianza assume una forma specifica per diverse distribuzioni.

Formule della Verosimiglianza:

1. Modello Normale: La verosimiglianza nel caso di una distribuzione normale è definita dalla densità di probabilità della distribuzione normale. Per una singola osservazione, la formula è: \[ L(y_i | \mu_i, \sigma^2) = \frac{1}{\sqrt{2\pi\sigma^2}} \exp\left(-\frac{(y_i - \mu_i)^2}{2\sigma^2}\right) \]

   Dove:

   • \(y_i\) è l’osservazione i-esima.
   • \(\mu_i\) è il valore atteso della variabile di risposta i-esima.
   • \(\sigma^2\) è la varianza.

2. Modello di Poisson: La verosimiglianza nel caso di una distribuzione di Poisson è definita come segue: \[ L(y_i | \lambda_i) = \frac{\lambda_i^{y_i} \exp(-\lambda_i)}{y_i!} \]

   Dove:

   • \(y_i\) è l’osservazione i-esima.
   • \(\lambda_i\) è il valore atteso della variabile di risposta i-esima.

Interpretazione: La devianza può essere interpretata come una misura della discrepanza tra il modello fitted e il modello null. Un valore di devianza più basso indica una migliore adattabilità del modello ai dati. Tuttavia, poiché la devianza è una misura assoluta, è spesso utilizzato il concetto di devianza residua, che è la devianza divisa per il numero di gradi di libertà del modello. Questo consente un confronto più equo tra modelli con differenti complessità.

Contributo delle Componenti: Nel contesto dei GLM, la devianza è spesso scomposta in tre componenti principali:

1. Modello Fitted Deviance (\(D_{\text{fitted}}\)): Misura la discrepanza tra il modello fitted e i dati osservati.

2. Null Deviance (\(D_{\text{null}}\)): Misura la discrepanza tra il modello null e i dati osservati.

3. Residual Deviance (\(D_{\text{residual}}\)): Rappresenta la devianza residua, cioè la discrepanza non spiegata dal modello fitted.

L’utilizzo di queste componenti permette di comprendere come la devianza è distribuita tra il modello fitted, il modello null e la devianza residua.

Confronto tra Modelli: Il test di devianza è spesso utilizzato per confrontare modelli alternativi. La differenza nella devianza tra due modelli segue approssimativamente una distribuzione chi-quadro sotto l’ipotesi nulla che i due modelli siano equivalenti. Questo test può essere utilizzato per valutare l’aggiunta di predittori al modello o per confrontare modelli con differenti specifiche di distribuzione della risposta.

Esempio di Devianza:

# Carichiamo un dataset di esempio in R
data(mtcars)

# Creiamo un modello di Poisson
model_poisson <- glm(vs ~ wt + hp, family = poisson, data = mtcars)

# Calcoliamo la devianza
deviance_value <- deviance(model_poisson)

cat("Devianza del Modello di Poisson:", deviance_value, "\n")

Devianza del Modello di Poisson: 9.969627 

In questo esempio, calcoliamo e visualizziamo la devianza residua di un modello di Poisson.

Residui

Funzione dei Residui nei GLM: Nei Modelli Lineari Generalizzati (GLM), i residui svolgono un ruolo cruciale nel valutare l’adattamento del modello ai dati e nel verificare la validità delle ipotesi dietro il modello. A differenza dei Modelli Lineari (LM), i residui nei GLM sono calcolati considerando la distribuzione della risposta specifica del modello.

1. Misurare la Bontà di Adattamento:
   • I residui nei GLM sono utilizzati per valutare la bontà di adattamento del modello. Se il modello si adatta bene ai dati, ci si aspetta che i residui abbiano una distribuzione che riflette la distribuzione della risposta specificata nel GLM.
2. Indicazioni sulla Struttura dei Dati:
   • Nei GLM, la scelta della distribuzione della risposta e della funzione di legame può variare in base alla natura dei dati. I residui forniscono indicazioni sulla struttura dei dati e sulla validità delle ipotesi del modello.
3. Diagnosticare Devianza:
   • La devianza, una misura della differenza tra il modello completo e uno più semplice, può essere diagnosticata attraverso i residui. Residui devianti e standardized deviance residuals sono spesso utilizzati per individuare modelli non adatti ai dati.

Differenze tra Residui nei GLM e nei LM: Le principali differenze tra i residui nei GLM e nei LM riguardano la distribuzione della risposta e la funzione di legame.

1. Distribuzione della Risposta:
   • Nei GLM, i residui sono calcolati tenendo conto della distribuzione della risposta specificata nel modello. Ad esempio, nei modelli di Poisson, i residui devono adattarsi alla distribuzione di Poisson.
   • Nei LM, i residui sono basati sull’assunzione che la risposta sia distribuita normalmente.
2. Funzione di Legame:
   • La funzione di legame nei GLM determina come il valore atteso della risposta è collegato alla combinazione lineare dei predittori. La scelta della funzione di legame influenza i residui.
   • Nei LM, la funzione di legame è identità, e i residui riflettono semplicemente la differenza tra i valori osservati e quelli previsti.

Esempio di Calcolo dei Residui in un Modello GLM:

suppressWarnings({
# Carichiamo un dataset di esempio in R con una distribuzione di Poisson
data(faithful, package = "datasets")

# Creiamo un modello di Poisson
modello_poisson <- glm(eruptions ~ waiting, family = poisson, data = faithful)

# Calcoliamo i residui devianti
residui_devianti <- residuals(modello_poisson, type = "deviance")

# Visualizziamo i primi 10 residui devianti
head(residui_devianti)
})

         1          2          3          4          5          6 
-0.2147663 -0.2970575 -0.1293332 -0.2561309 -0.0479040  0.3822210 

Analisi dei Residui: I residui devianti riflettono le differenze tra i valori osservati e quelli previsti in termini della devianza del modello. Un residuo deviante elevato indica che l’osservazione contribuisce in modo significativo alla devianza complessiva del modello, indicando un’eventuale influenza o deviazione dal modello.

• Residui positivi indicano che l’osservazione ha contribuito più del previsto alla devianza.
• Residui negativi indicano che l’osservazione ha contribuito meno del previsto alla devianza.
• Residui pari a zero indicano una perfetta adattabilità dell’osservazione al modello.

Esempio

Per questo esempio, useremo il dataset di esempio “mtcars” di R per creare un modello di Poisson utilizzando un modello generalizzato lineare (GLM). Lo scopo del modello sarà prevedere il numero di cilindri (“cyl”) in base alle altre variabili presenti nel dataset.

# Caricamento del dataset "mtcars"
data(mtcars)

# Esploriamo le prime righe del dataset
head(mtcars)

# Creiamo un modello di Poisson per prevedere il numero di cilindri in base alle altre variabili
modello_glm <- glm(cyl ~ mpg + disp + hp + drat + wt + qsec + vs + am + gear + carb, 
                   data = mtcars, family = poisson)

# Visualizziamo il summary del modello
summary(modello_glm)

Call:
glm(formula = cyl ~ mpg + disp + hp + drat + wt + qsec + vs + 
    am + gear + carb, family = poisson, data = mtcars)

Deviance Residuals: 
     Min        1Q    Median        3Q       Max  
-0.42777  -0.16289   0.01851   0.15996   0.40340  

Coefficients:
              Estimate Std. Error z value Pr(>|z|)
(Intercept)  3.1253653  2.3664672   1.321    0.187
mpg         -0.0058390  0.0351590  -0.166    0.868
disp         0.0006115  0.0025727   0.238    0.812
hp           0.0002488  0.0033108   0.075    0.940
drat        -0.0981988  0.2466800  -0.398    0.691
wt          -0.0366284  0.3070279  -0.119    0.905
qsec        -0.0406992  0.1219151  -0.334    0.739
vs          -0.1092921  0.3282386  -0.333    0.739
am          -0.1126824  0.3368595  -0.335    0.738
gear        -0.0517634  0.2304110  -0.225    0.822
carb         0.0280157  0.1270379   0.221    0.825

(Dispersion parameter for poisson family taken to be 1)

    Null deviance: 16.574  on 31  degrees of freedom
Residual deviance:  1.310  on 21  degrees of freedom
AIC: 139.97

Number of Fisher Scoring iterations: 4

# Analizziamo le variabili indipendenti
par(mfrow = c(2, 2))
plot(modello_glm)

# Eseguiamo l'analisi della varianza (ANOVA)
anova_result <- anova(modello_glm, test = "Chi")

# Visualizziamo la tabella ANOVA
print(anova_result)

Analysis of Deviance Table

Model: poisson, link: log

Response: cyl

Terms added sequentially (first to last)

     Df Deviance Resid. Df Resid. Dev  Pr(>Chi)    
NULL                    31    16.5743              
mpg   1  12.2921        30     4.2822 0.0004549 ***
disp  1   1.1779        29     3.1043 0.2777770    
hp    1   0.3231        28     2.7812 0.5697385    
drat  1   0.3586        27     2.4225 0.5492712    
wt    1   0.2091        26     2.2135 0.6475086    
qsec  1   0.4558        25     1.7577 0.4996039    
vs    1   0.1329        24     1.6248 0.7154335    
am    1   0.2440        23     1.3808 0.6213230    
gear  1   0.0221        22     1.3587 0.8817565    
carb  1   0.0487        21     1.3100 0.8253848    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

# Calcoliamo l'R-squared del modello
r_squared <- 1 - (modello_glm$deviance / modello_glm$null.deviance)
cat("R-squared:", r_squared, "\n")

R-squared: 0.9209632 

# Effettuiamo previsioni su nuovi dati (per esempio, le prime 5 osservazioni del dataset)
nuovi_dati <- mtcars[1:5, ]
previsioni <- predict(modello_glm, newdata = nuovi_dati, type = "response")
cat("Previsioni per le prime 5 osservazioni:\n", previsioni, "\n")

Previsioni per le prime 5 osservazioni:
 5.879263 5.693355 4.305739 5.683563 7.787923 

In questo esempio, abbiamo creato un modello di Poisson utilizzando il numero di cilindri come variabile dipendente e le altre variabili del dataset “mtcars” come variabili indipendenti. Successivamente, abbiamo eseguito un’analisi completa del modello, compresi il summary, la visualizzazione delle variabili indipendenti, l’analisi della varianza (ANOVA), il calcolo dell’R-squared e la previsione su nuovi dati.

[[Torna all’ Indice]]